% GAUTIER LE BIHAN - 2020
% Replication files for �Shocks vs Menu Costs: Patterns of Price Rigidity in an Estimated Multi-Sector Menu-Cost Model� �Review of Economics and Statistics
%
% This code Plots  Appendix figure H based on stored  simulated data

load ..\..\ident_canova_sala_productivity\mat_results_CS ;

param=mat_results_CS;

moments=[ 0.0837  0.6737  0.0438  0.0698  4.4616];
 
actual_std=[ 0.0008  0.0042  0.0005  0.0010  0.1057]
actual_var= actual_std.* actual_std*1000;

 
m_target = moments';
m_scale=[actual_var'];
s=size(param,1);

for jj=1:s;
m_model=param(jj,6:10)';
G = m_target - m_model;
% % Weighting matrix W
size_weights = max(size(m_target));
scale = (ones(size_weights,1).*[4;1;1;1;1])./m_scale; 
% % take probability of default with 'double importance'
 W = eye(size_weights).*(ones(size_weights,1)*scale');
 clear scale size_weights
% 
 dist_sqrd = G.*G;
 W_vector = [W(1,1);W(2,2);W(3,3);W(4,4);W(5,5)];
dist_weighted = dist_sqrd.*W_vector;
% % Calculate new SMM estimator:
% 
 x(jj,1) = G'*W*G;
end;
param=[param x];

dist=param(:,11);





mylevel=-1*[ 0.004 0.013 0.033  0.05 0.07499]% 0.099]
%mylevel=-1*[ 0.008 0.012  0.035  0.05 0.07]
%mylevel=-1*[ 0.008 0.015  0.035  0.099]

%param1=param((dist<0.1 & dist>0 ) ,:)
param1=param((dist<0.065 & dist>0 ) ,:);

ncontour=4
paramopt=param((dist==min(dist((dist>0),1)) ),:);


param2=param1((param1(:,1)~= paramopt(1)),:);
dist2=param2(:,11);
paramopt2=param2((dist2==min(dist2) ),:);
%param0=[0.0501    0.0373   0.0381   0.6946];

M=[param1(:,1) param1(:,2) -1*param1(:,11)];
x = linspace(min(M(:,1)), max(M(:,1)), size(M,1)+1);
y = linspace(min(M(:,2)), max(M(:,2)), size(M,1)+1);
[X,Y] = meshgrid(x, y);
Z = griddata(M(:,1), M(:,2), M(:,3), X, Y);
figure(1)
subplot(2,3,1)
[C,h]=contourf(X, Y, Z, ncontour);
h.LevelList=mylevel;

map =  [
%ext�rieur   
 0.8 0.8 0.8
%2nd stage
      0.6 0.6 0.6
      0.6 0.6 0.6
      %3rd stage
      0.4 0.4 0.4 
      0.4 0.4 0.4 
      %Center
      0 0 0
       %0 0 1 
    
    ];
colormap(map)
%colormap(pink(5))
cmap = colormap
x0=paramopt(:,1);
y0=paramopt(:,2);
%return;
c=[1 0 0]
sz=50;
hold on
scatter(x0, y0,sz,'r', 'filled');
hold off
hold on
scatter(paramopt2(:,1), paramopt2(:,2),sz,c,'g', 'filled');
hold off
xlabel('\lambda') % x-axis label
ylabel('\mu') % y-axis label

%contour(X, Y, Z, 5, 'ShowText','on')
%axis equal
%grid

%figure(2)
%surf(X, Y, -1*Z)
%contour(X, Y, Z, 5, 'ShowText','on')
%axis equal
%grid

M=[param1(:,1) param1(:,3) -1*param1(:,11)];
x = linspace(min(M(:,1)), max(M(:,1)), size(M,1)+1);
y = linspace(min(M(:,2)), max(M(:,2)), size(M,1)+1);
[X,Y] = meshgrid(x, y);
Z = griddata(M(:,1), M(:,2), M(:,3), X, Y);
subplot(2,3,2)
[C,h]=contourf(X, Y, Z, ncontour);
h.LevelList=mylevel;
%colormap(pink)
map =  [
%ext�rieur   
 0.8 0.8 0.8
%2nd stage
      0.6 0.6 0.6
      0.6 0.6 0.6
      %3rd stage
      0.4 0.4 0.4 
      0.4 0.4 0.4 
      %Center
      0 0 0
       %0 0 1 
    
    ];
% map =  [
% %ext�rieur   
% 0 0.0 0.0
% %2nd stage
%       1 0 0.
%       1 0 0.
%       %3rd stage
%       00. 1 0.
%       00. 1 0.
%       %Center
%       00. 0 1
%        %0 0 1 
%     
%     ];
%surf(peaks)
colormap(map)
%colormap(pink(5))
x0=paramopt(:,1);
y0=paramopt(:,3);
%x1=.0501;
%y1=.0373;
sz=60
c=[1 0 0]
hold on
scatter(x0, y0,sz,'r', 'filled');
hold off
hold on
scatter(paramopt2(:,1), paramopt2(:,3),sz,c,'g', 'filled');
hold off
xlabel('\lambda') % x-axis label
ylabel('\sigma') % y-axis label


%contour(X, Y, Z, 5, 'ShowText','on')
%axis equal
%grid

%figure(4)
%surf(Y, X, -1*Z)


M=[param1(:,1) param1(:,4) -1*param1(:,11)];
x = linspace(min(M(:,1)), max(M(:,1)), size(M,1)+1);
y = linspace(min(M(:,2)), max(M(:,2)), size(M,1)+1);
[X,Y] = meshgrid(x, y);
Z = griddata(M(:,1), M(:,2), M(:,3), X, Y);
subplot(2,3,3)
%contour(X, Y, Z, ncontour)
[C,h]=contourf(X, Y, Z, ncontour);
h.LevelList=mylevel;
map =  [
%ext�rieur   
 0.8 0.8 0.8
%2nd stage
      0.6 0.6 0.6
      0.6 0.6 0.6
      %3rd stage
      0.4 0.4 0.4 
      0.4 0.4 0.4 
      %Center
      0 0 0
       %0 0 1 
    
    ];
% map =  [
% %ext�rieur   
% 0 0.0 0.0
% %2nd stage
%       1 0 0.
%       1 0 0.
%       %3rd stage
%       00. 1 0.
%       00. 1 0.
%       %Center
%       00. 0 1
%        %0 0 1 
%     
%     ];
%surf(peaks)
colormap(map)
%colormap(pink(5))
%colormap(pink)
x0=paramopt(:,1);
y0=paramopt(:,4);
%x1=.0501;
%y1=.0381;
sz=60
c=[1 0 0]

hold on
scatter(x0, y0,sz,c, 'filled');
hold off
hold on
scatter(paramopt2(:,1), paramopt2(:,4),sz,c,'g', 'filled');
hold off
xlam=[0.046 0.054]
ysig=[0.032 0.038]
%axis([xlam ysig]);
%hold on
%scatter(x1, y1,sz,c, 'filled');
%hold off
xlabel('\lambda') % x-axis label
ylabel('\rho') % y-axis label

%figure(6)
%surf(Y, X, -1*Z)




M=[param1(:,2) param1(:,3) -1*param1(:,11)];
Z = griddata(M(:,1), M(:,2), M(:,3), X, Y);
x = linspace(min(M(:,1)), max(M(:,1)), size(M,1)+1);
y = linspace(min(M(:,2)), max(M(:,2)), size(M,1)+1);
[X,Y] = meshgrid(x, y);
Z = griddata(M(:,1), M(:,2), M(:,3), X, Y);
subplot(2,3,4)
%contour(X, Y, Z, ncontour)
[C,h]=contourf(X, Y, Z, ncontour);
h.LevelList=mylevel;
map =  [
%ext�rieur   
 0.8 0.8 0.8
%2nd stage
      0.6 0.6 0.6
      0.6 0.6 0.6
      %3rd stage
      0.4 0.4 0.4 
      0.4 0.4 0.4 
      %Center
      0 0 0
       %0 0 1 
    
    ];
% map =  [
% %ext�rieur   
% 0 0.0 0.0
% %2nd stage
%       1 0 0.
%       1 0 0.
%       %3rd stage
%       00. 1 0.
%       00. 1 0.
%       %Center
%       00. 0 1
%        %0 0 1 
%     
%     ];
%surf(peaks)
colormap(map)
%colormap(pink(5))
%colormap(pink)
x0=paramopt(:,2);
y0=paramopt(:,3);
%x1=.0501;
%y1=.0381;
sz=60
c=[1 0 0]

hold on
scatter(x0, y0,sz,c, 'filled');
hold off
hold on
scatter(paramopt2(:,2), paramopt2(:,3),sz,c,'g', 'filled');
hold off
xlabel('\mu') % x-axis label
ylabel('\sigma') % y-axis label

%figure(8)
%surf(Y, X, -1*Z)

M=[param1(:,2) param1(:,4) -1*param1(:,11)];
x = linspace(min(M(:,1)), max(M(:,1)), size(M,1)+1);
y = linspace(min(M(:,2)), max(M(:,2)), size(M,1)+1);
[X,Y] = meshgrid(x, y);
Z = griddata(M(:,1), M(:,2), M(:,3), X, Y);
subplot(2,3,5)
%contour(X, Y, Z, ncontour)
[C,h]=contourf(X, Y, Z, ncontour);
h.LevelList=mylevel;
map =  [
%ext�rieur   
 0.8 0.8 0.8
%2nd stage
      0.6 0.6 0.6
      0.6 0.6 0.6
      %3rd stage
      0.4 0.4 0.4 
      0.4 0.4 0.4 
      %Center
      0 0 0
       %0 0 1 
    
    ];
% map =  [
% %ext�rieur   
% 0 0.0 0.0
% %2nd stage
%       1 0 0.
%       1 0 0.
%       %3rd stage
%       00. 1 0.
%       00. 1 0.
%       %Center
%       00. 0 1
%        %0 0 1 
%     
%     ];
%surf(peaks)
colormap(map)
%colormap(pink(5))
%colormap(pink)
x0=paramopt(:,2);
y0=paramopt(:,4);
%x1=.0501;
%y1=.0381;
sz=60
c=[1 0 0]

hold on
scatter(x0, y0,sz,c, 'filled');
hold off
hold on
scatter(paramopt2(:,2), paramopt2(:,4),sz,c,'g', 'filled');
hold off
xlabel('\mu') % x-axis label
ylabel('\rho') % y-axis label

%figure(10)
%surf(Y, X, -1*Z)

M=[param1(:,3) param1(:,4) -1*param1(:,11)];
x = linspace(min(M(:,1)), max(M(:,1)), size(M,1)+1);
y = linspace(min(M(:,2)), max(M(:,2)), size(M,1)+1);
[X,Y] = meshgrid(x, y);
Z = griddata(M(:,1), M(:,2), M(:,3), X, Y);
subplot(2,3,6)
%contour(X, Y, Z, ncontour)
[C,h]=contourf(X, Y, Z, ncontour);
h.LevelList=mylevel;
map =  [
%ext�rieur   
 0.8 0.8 0.8
%2nd stage
      0.6 0.6 0.6
      0.6 0.6 0.6
      %3rd stage
      0.4 0.4 0.4 
      0.4 0.4 0.4 
      %Center
      0 0 0
       %0 0 1 
    
    ];
% map =  [
% %ext�rieur   
% 0 0.0 0.0
% %2nd stage
%       1 0 0.
%       1 0 0.
%       %3rd stage
%       00. 1 0.
%       00. 1 0.
%       %Center
%       00. 0 1
%        %0 0 1 
%     
%     ];
%surf(peaks)
colormap(map)
%colormap(pink(5))
%colormap(pink)
x0=paramopt(:,3);
y0=paramopt(:,4);
%x1=.0501;
%y1=.0381;
sz=60
c=[1 0 0]


hold on
scatter(x0, y0,sz,c, 'filled');
hold off
hold on
scatter(paramopt2(:,3), paramopt2(:,4),sz,c, 'g','filled');
hold off
xlabel('\sigma') % x-axis label
ylabel('\rho') % y-axis label
 print('..\..\..\figures\Canova_rho.pdf','-dpdf','-fillpage')

%figure(12)
%surf(Y, X, -1*Z)










figure(13);
%P0
subplot(2,3,1);
x=param1(:,1);
y=param1(:,2);
x0=paramopt(:,1);
y0=paramopt(:,2);


scatter(y, x);
hold on
scatter(y0, x0);
hold off
x1=[0.01  0.06]
y1=[0.04  0.06]


axis([x1 y1])
ylabel('\mu') % x-axis label
xlabel('\lambda') % y-axis label

subplot(2,3,2);
x=param1(:,1);
y=param1(:,3);
x0=paramopt(:,1);
y0=paramopt(:,3);


scatter(y, x);
hold on
scatter(y0, x0);
hold off
x1=[0.025  0.06]
y1=[0.04  0.06]


axis([x1 y1])
ylabel('\lambda') % x-axis label
xlabel('\sigma') % y-axis label


subplot(2,3,3);
x=param1(:,2);
y=param1(:,3);
x0=paramopt(:,2);
y0=paramopt(:,3);


scatter(y, x);
hold on
scatter(y0, x0);
hold off
x1=[0.02  0.06]
y1=[0.02  0.045]


axis([x1 y1])
ylabel('\mu') % x-axis label
xlabel('\sigma') % y-axis label


subplot(2,3,4);
x=param1(:,1);
y=param1(:,4);
x0=paramopt(:,1);
y0=paramopt(:,4);


scatter(y, x);
hold on
scatter(y0, x0);
hold off
x1=[0.5  0.8]
y1=[0.04  0.06]


axis([x1 y1])
ylabel('\mu') % x-axis label
xlabel('\rho') % y-axis label


subplot(2,3,5);
x=param1(:,2);
y=param1(:,4);
x0=paramopt(:,2);
y0=paramopt(:,4);


scatter(y, x);
hold on
scatter(y0, x0);
hold off
x1=[0.5  0.8]
y1=[0.02  0.045]


axis([x1 y1])
ylabel('\mu') % x-axis label
xlabel('\rho') % y-axis label


subplot(2,3,6);
x=param1(:,3);
y=param1(:,4);
x0=paramopt(:,3);
y0=paramopt(:,4);


scatter(y, x);
hold on
scatter(y0, x0);
hold off
x1=[0.5  0.8]
y1=[0.02  0.06]


axis([x1 y1])
ylabel('\sigma') % x-axis label
xlabel('\rho') % y-axis label

 print('..\..\..\figures\Canova_ident_rho.pdf','-dpdf')
 
